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^ ■ Abstract 

^D I The density profiles and other quantities of physical interest for spherically 

^^ ' symmetric systems are computed by assuming that a collisionless stellar gas 

^^ 
Mh| may relax to the non-Gaussian power law distribution suggested by the nonex- 

o, 

H , tensive kinetic theory. There are two different classes of solutions. The first 

' class behaves like a subset of the polytropic Lane-Emden spheres, whereas 

k>>( , the second one corresponds to a transition between two different polytropic 

H , 

C^ I indices. Unlike the isothermal Maxwellian sphere, the total mass and sizes of 

both classes are finite for a large range of the nonextensive q-parameter. 
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I. INTRODUCTION 

It is widely assumed that self-gravitating stellar systems like globular clusters are com- 
pletely or almost completely relaxed because they are changing very slowly, or equivalently, 
the characteristic evolution time scale is too long. In particular, this means that stellar 
encounters combined with the action of some long range kinetic processes like "phase mix- 
ing" and "violent relaxation" have driven their velocity distributions toward a collisionless 
kinetic steady state (Spitzer & Harm 1958, King 1962, Michie 1961, 1963, Lynden-Bell 1967, 
Spitzer 1987). Although considering the numerous and important contributions to this field, 
an overlook in the recent literature suggests that the equilibrium (or near-equilibrium) distri- 
bution of these relaxed systems is still not firmly established (beyond doubt), and basically 
remains as an open question (Drukier et al. 1992, Meylan & Heggie 1997, Katz & Okamoto 
2000, Chavanis 2002). 

The phase-space density for spherically symmetric systems was initially described by the 
equilibrium Maxwell-Boltzmann (MB) distribution 

/ ^x-s/s / ^v'^ + $(r)\ , , 

f{r,v) = pi{2na') exp ( - ^ ^^ ^ M (1) 

where $ is the potential per unit mass, pi is the density in the absence of the force field 
and a is the velocity dispersion. However, the computation of the associated mass density 
gave rise to a very serious problem, namely: the total mass is infinite (Chandrasekhar 1960, 
Binney & Tremaine 1994). This undesirable feature is more easily seen from the particular 
solution, p(r) = o"^/27rG'r^ (singular isothermal sphere) which predicts an infinite total mass. 
Many ad hoc phenomenological distributions were proposed to solve this problem with basis 
on Jeans' theorem: any steady state solution of the collisionless Boltzmann (or Vlasov) 
equation depends on the phase space coordinates only through isolating integrals of motion 
(Jeans 1929). 

A different route was initiated by Spitzer and Harm (1958), and further worked out 
by Michie (1961, 1963) and King (1962, 1965, 1966). It was recognized that star clusters 



cannot be exactly described by the expected MB equilibrium distribution. However, it 
was also assumed that relaxation processes always drive the stellar distribution as far as 
it can to a quasi-MB final state which should be determined, for instance, by solving the 
non-equilibrium Fokker-Planck equation. This line of inquiry lead to some interesting and 
consistent results among them the lowered isothermal sphere, as well as the anisotropic 
models proposed by Michie. 

In this article we consider a different approach. With basis on the q-power law equilibrium 
distribution we determine the radial and projected density profiles for two large classes of 
isothermal stellar systems. It should be recalled that the nonextensive treatment for a stellar 
collisionless systems was first considered by Plastino & Plastino (1993) through a variational 
principle where the Tsallis (1988) entropy formula was maximized taking into account the 
constraints imposed by the total mass and energy density. Here, we consider directly the 
nonextensive distribution which follows naturally from the kinetic equilibrium q-entropy 
formula (Lima et al. 2001). An attractive feature of this power law distribution is that the 
models are analytically tractable in such a way that a detailed comparison with the standard 
Maxwell-Boltzmann approach is immediate. As we shall see, there are at least two classes 
of solutions with finite mass and radius. Actually, the first class is not new, however we 
discuss some related subtleties not considered in the quoted articles. 

The paper is structured as follows. Next section we set up the basic equations and discuss 
the density profiles to the first class of nonextensive spherically symmetric self-gravitating 
stellar systems based on Tsallis' distribution. Section 3 discuss the main features of the 
trunctated models, and, in section 4, we resume briefiy the main results. 

II. POWER LAW STELLAR SPHERES 



Let us now consider the following non-Gaussian velocity distribution 
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where the g-parameter quantifies to what extent the distribution departs from the standard 
Gaussian form (Silva et ah 1998, Lima at al. 2001, Mendes & Tsalhs 2001, Kaniadakis 
2001). Actually, it is a power law if g 7^ 1, whereas for g = 1 it reduces to the Maxwellian 
function. This result follows directly from the known identity limd^o(l + dy)d = exp(y) 
(Abramowitz & Stegun 1972). The quantity Cq is a g-dependent normalization constant 
given by (Lima et al. 2002) 

C, = il-qy^^'-^)i^^)^^ forg<l (3) 

and 

C, = il- If' (^) T^^W^ /- 1 ^ 1' (4) 

[l-q 2) 

which reduce to the expected result in the limit q = 1. For instance, for g < 1 we define 
z = (1— g)^^ so that 2; ^ 00 if g goes to unity. In addition, from the identity lim|z|^ooZ~'^r(a+ 
z)/r(z) = 1 (Abramowitz & Stegun 1972), one may see that Ci = 1, thereby showing that 
the normalization of the standard Maxwellian distribution is recovered in this limit. 

The power law distribution (2) can be deduced from two simple requirements: (i) 
isotropy of the velocity space, and (ii) a suitable nonextensive generalization of the Maxwell 
factorizability condition, or equivalently, the assumption that in this enlarged framework 
/('^) 7^ f {vx) f {vy) f {Vz) ■ It was also shown that for g > 0, the above distribution func- 
tion satisfies a generalized H-theorem, and its reverse has also been proved, that is, the 
collisional nonextensive equilibrium is given by the Tsallis' power law velocity distribution 
(Lima et al. 2001). In the last few years, several applications of this equilibrium power 
law velocity distribution have been done in many disparate branches of physics (Liu et al. 
1994, Bhogosian 1996, Lima et al. 2000, Tsallis et al. 2001). In the astrophysical context, 
Taruya and Sakagami (2002), investigated in detail the problem of instability (gravothermal 
catastrophe) working in the so-called (u,v) plane (Milne's variables). Further, they shown 
that the king model fails to match the simulated distribution function, especially at cut-off 
scales (Taruya & Sakagami 2003). Even the Jeans gravitational instability criterion for a 



collisionless system was recently discussed in this enlarged framework (Lima et al 2002). In 
particular, for power law distribution with cut-off, it was shown that such a system presents 
instability even for wavelengths of the disturbance smaller than the standard critical Jeans 
value. For a selfgravitating system, the nonextensive phase-space density reads (Lima et al. 

2002a) 
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which can also be obtained by integrating the Vlasov equation. As should be expected, the 
above expression reduces to equation (1) in the limiting case q = 1- 

Let us now consider the velocity distribution (5) to build a set of stellar systems whose 
structure is sustained by their respective gravitational field. Following standard lines, the 
granularity of the star system is ignored, the gravitational potential is assumed to be a slowly 
varying function of position, and any change in the physical properties due to collisions or 
evolution of the stars are neglected. The radial dependence of the mass density is obtained 
by integrating (5) over all allowed velocities 
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whereas the gravitational potential must be determined from Poisson's equation 

}_d_ 2d4_. ^ ^^Q 
r^ dr dr 



(6) 
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In order to simplify further the integral (6), it is convenient to introduce the dimensionless 
energy per unit mass of a star 






(8) 



and from the argument of the power law in (5), this quantity is restricted by e < j^, in 
order to guarantee a definite real valued distribution function (unless explicitly stated, in 
what follows we consider only the case q < 1). Now, inserting e into (6) one has 

4:npiCg fih; 
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and a simple integration furnishes 
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p = pjl-(l-g)^]^(^ (10) 

where the density scale is expressed in terms of the complete Beta function by Pq = 
pi(l — g)~2_B(|, 1^). In the limit q = I one finds Pq = pi with the power-law becom- 
ing the exponential, and as expected the standard Maxwellian result is readily recovered. 
Substituting (10) into Poisson's equation, and introducing the pair of dimensionless quanti- 
ties defined by a; = r/vo and = —a^'^(f){r), where r^ = y^irGp/a'^, we find the Lane-Emden 
type equation 

whose importance for stellar dynamics is largely known. First, we notice that in the limit 
q = I the above equation reduces to 

-M^'^) = -e' (12) 

which is exactly the Lane-Emden form for a Maxwellian isothermal sphere (Chandrasekhar 
1960, Spitzer 1987, Binney & Tremaine 1994). The solution of this equation subjected 
to the boundary conditions 9 = dO/dx = at a; = has been numerically computed by 
Chandrasekhar and coworkers (see Chandrasekhar 1960 and references there in). Some 
important theorems has been proven both for the Maxwellian isothermal sphere and the 
Lane-Emden equation. Such a study can also be extended to equation (11) adopting the 
same boundary conditions. 

Now, it should be recalled that the gravitational potential is defined up to a constant 
value in an arbitrary level. Therefore, one may also choose 0o in such a way that (f) -^ (p — (po 
subjected to the restriction that both the potential and the density goes to zero at infinite. 
With that choice, the arbitrary constant for q ^ 1 is (po = 1/(1 ^5) ^^^ from (10) we obtain 
the following expression 

P = p,e^^ (13) 



whereas the differential equation (11) becomes 
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which is the canonical form of a polytropic Lane-Emden equation of index 

In figure 1 we display the profile p/pq for some values of the q parameter. The labels 
under each curve correspond to the value lOg and since the density profiles are identical to 
the classical polytropic solution we obtain the known restriction that the total mass diverges 
above the polytropic index n = 5 since the Lane-Emden spheres have finite mass only for 
n < 5. In terms of the g-parameter this condition is simply translated as g < 5/7. On the 
other hand, from the second law of thermodynamics we know that q > (Lima et al 2001). 
Therefore, models with finite mass and sizes are obtained if the nonextensive parameters lies 
on the interval < g < 5/7. Such a limit was previously obtained by Plastino and Plastino 
(1993) maximizing Tsallis entropy. The unique difference is that it was referred to as 9/7 
because the distribution with cut-off was written as q-1 instead of 1-q. In fact, substituting 
q ^ 2 — q we see that 5/7 goes to 9/7, with the latter becoming a lower bound as should 
be expected. Note that the physical lower bound constraint (g > 0) means that n > 5/2. 
Therefore, Tsallis distributions whose total energy has a cut-off ranges only a half of all 
possible Lane-Emden distributions with finite mass. However, since the above Lane-Emden 
type equations are still valid for g > 1, one may see that n > means g > 5/3, and in the 
limit g — ^ oo the index n -^ 3/2. Hence, we conclude that Tsallis power law statistics do not 
describe the class of polytropic distributions contained on the range 3/2 < n < 5/2 since 
this range requires negative values of g. These results can be seen directly in the plane (n,q) 
represented in figure 2. 

At light of such results, one may be tempted to conclude that a subset of Lane-Emden 
stellar polytropes can be regarded to as nonextensive collisionless isothermal spheres, thereby 
giving a kinetic justification to this class of stellar configurations. It is argued here, however, 
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that such an identification is misleading since polytropic models are defined by one of their 
equation of state, namely: pp^'^ = const, or equivalently, Tp^^'^ = const, where 7 = 1 + 1/n is 
the polytropic index. Hence, the temperature in polytropic models depends on the position, 
whereas the class of nonextensive spheres (including the Maxwellian case) has constant 



temperature which can formally be defined by the velocity dispersion a = JkBT/m. 

In figure 3 we display the projected mass density S for the same set of the q-parameter 
illustrated in figure 1. These profiles illustrates the surface density profile that should be 
obtained in an astronomical object described by this class of Tsallis models. 

It is worth notice that Tsallis isothermal spheres have a fraction of the stars with positive 
energy and even so they are bounded to the main structure whose total mass is finite for 
q < 5/7. In this concern, one may ask what happens if the energy is constrained by £ < 
as usually adopted for bounded structures. Physically, one would expect finite structures 
for a larger range of the q-parameter. This problem it will be discussed next section. 

III. TRUNCATED NONEXTENSIVE STELLAR MODELS 

Let us now consider a class of nonextensive stellar models which is obtained by assuming 
that only objects with e < are present in the distribution. It can be viewed as a sim- 
ple generalization of the truncated model proposed by Wooley (1954). The corresponding 
density profile as a function of the gravitational potential assumes the form (Cf. equation 
(9)) 



p = ^^ Ji-a-^M^^is-^M (16) 



where integration limits are defined by equation (8). As one may check, this integral can be 
expressed in terms of hypergeometric functions in the form 

^1 3 ?7* 2 1 - g 2 r]* 

where rf^Lp) = 1 — (1 — q)Lp, r]^ = (^Z) ? ^^(^ '■P = -^- In order to work out a numerical 
solution it is more convenient to express this density profile in a form where the integration 



limits are fixed as 

p = C,{-vfl^ (\l - (1 - g)v^x]^(l - xY^dx (18) 

Jo 

and one may verify tliat at tlie outer lialo region, wliere \ip\ « 1, tlie density profile behaves 
like 

p^p,(-^)3/2 k|«l (19) 

corresponding to a polytropic spherical profile having Lane-Emden index n = 3/2. In the 
central region, where we might adopt the approximation \ip\ >> 1, the density profile reduces 
to the expression 

p = p^{-^)W^ (20) 

which corresponds to a Lane-Emden index n = (5 — 3g)/2(l — q). Therefore, if all the 
stars has velocity smaller than the escape velocity, the general density profile in terms of 
the gravitational potential can be described as a smooth transition between two different 
polytropic spheres, where the case n = 3/2 behaves like an attractor for the outer halo region. 
This behaviour is illustrated in figure 4. As a further consequence all these truncated models 
have finite mass extent due to the very steep density profile dominating the external region. 
The density profile of these truncated models can be easily obtained by a numerical 
integration of the the Poisson equation subject to the auxiliary condition given by equation 
17. The integration begins at r = by specifying a value ipo for the gravitational potential 
at the center of the mass distribution of the model. Knowing the potential we estimate 
the central density predicted by equation 17. Due to the continuity equation the density 
gradient at the center is null, and the same happens with the gradient of the gravitational 
potential. Once the density and the gradient of the potential are determined we apply a 
discretization algorithm to the second order Poisson equation for a grid of radial points. The 
algebraic approximation is used to predict the gravitational potential, and its gradient, in 
the next point of the grid. The whole process is iterated by solving the density in the next 
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point of the grid and applying recursively the discretization of the Poisson equation. The 
overall numerical integration stop at the point where the gravitational potential is null. At 
that point the corresponding density is also null and we reach the surface of the truncated 
model. 

In figure 5 we present the result of applying this integration algorithm for the cases q = 1, 
0.9 and 0.8. In each panel the continuous line corresponds to the non truncated model. 
This model corresponds to the asymptotic case when the central gravitational potential is 
sufficiently high to sample the polytropic structure. In the particular case of g = 1 we 
recover the classical isothermal sphere. On the other extreme we have the situation where 
the central gravitational potential tends to zero. In that case we can see from figure 4 that 
the structure is described by a n = 3/2 politrope independently from the q value. For 
that reason the density profiles for the three models in figure 5 are exactly the same when 
ipo -^ 0. For intermediate values of the central potential the models tend to be quite different 
depending on the exact value of q. In that figure we have chosen in each panel two values 
of ipo corresponding to intermediate density profiles. 

The density profiles can be integrated along the line of sight so that we can estimate 
the projected mass density. Assuming that the light trace the mass and that the mass to 
luminosity ratio is constant we can scale this to the surface brightness. In figure 6 we present 
this surface brightness as a function of the r^/^ radial coordinate. The symbols are the same 
as in the previous figure. Again we can see that large values of ipo tend to be closer to 
the nontruncated polytropic models as indicated by the continuous line. When ipo ^ the 
models are described by the surface mass density obtained from the n = 3/2 polytrope. 

IV. FINAL COMMENTS 

Two simple applications of Tsallis' power law kinetic distribution have been discussed. 
The proposed kinetic models represent the mass distribution of coUisionless gravitational 
stellar systems. The use of the Tsallis distribution introduces an extra parameter measuring 
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to what extent the velocity distribution departs from the standard Maxwell-Boltzmann law 
which is normally used to represent these objects. It has been proposed as a viable and well 
grounded alternative to the standard equilibrium approach in the presence of long range 
forces as happens in the astrophysical context. 

As we have seen, the Poisson equation combined with Tsallis distribution usually repro- 
duces spherically symmetric structures resembling the classical polytropic spheres. If the 
natural cut-off is imposed {q < 1), the polytropic index n is closely related to this Tsallis 
parameter. The basic result is that models with q < 5/7 have finite mass while above this 
limit the mass is divergent. The radial and projected density profiles were obtained by 
solving numerically the Poisson equation for a large range of the nonextensive parameter. 
Moreover, if objects with positive energies are excluded, one can build a larger set of trun- 
cated models with finite mass and sizes. For this class of models the corresponding profiles 
were also numerically determined (see figures 5 and 6). For both cases, the predicted profiles 
should be compared with that ones observed for globular clusters. In particular, since energy 
truncated models are presumed to describe the relaxed state of clusters, it is interesting to 
investigate their connection with possible nonextensive extensions of the spherically sym- 
metric Michie-King models. Finally, we also remark that the possibility to include elliptical 
galaxies in this generalized framework cannot be discarded. These issues will be discussed 
in a forthcoming communication. 
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FIG. 1. Spatial density profile for different values of q. For q < 5/7 all profiles have finite 
extent and finite mass, while for 5/7 < q < 1 the models extend to infinite and the total mass 
diverges. The limiting value q = 5/7 corresponds to n = 5 polytropic model. 
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FIG. 2. The q index of the Tsallis isothermal description versus the polytropic index n. In the 
left region we have < g < 1 and the corresponding polytropic models are restricted to the interval 
5/3 < n < oo. The polytropic models < n < 3/2 are described by Tsallis models having q > 5/3. 
Note that polytropic models defined by 3/2 < n < 5/2 do not have a corresponding Tsallis model. 
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FIG. 3. Projected mass distribution for different values of q, obtained by numerical integration 
of the profiles presented in figure 1. 
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FIG. 4. A class of truncated models can be obtained by extracting from the density distribution 
all objects with positive energy. In that case the density profile changes and can be described by 
a superposition of two limiting polytropic models. 
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FIG. 5. Density profiles predicted for the truncated model corresponding to q' = 1, 0.9 and 
0.8. The continuous lines in each panel describes the nontruncated polytrope corresponding to the 
limiting case of large c^o- When ipo ^ the models approximate the limiting n = 3/2 polytrope 
independently of the value of q. 







q = 1 .0 




\\, 


____ %^-s 




■ %^-'^ 


\ ■■- ^v^^ 


- - (|>|., -> 


- : ^ \^^^^ 


- 






q = 0.9 




' ^\ ^x 


(|>„ = -15 


■ ■■ ^v 


% = -"^ 


' '■-. ^-^ 


- - % -> 


- ' . "'^-x^^^^ 


- 






q = 0.8 




\ -'\ ^\ 


___. % = -30 


■ \ ^\ 


(|>(, = -10 


^v 


- - <|>o -> 


- \ \\^\^ 


- 





FIG. 6. Projected mass density for each of the density profiles of figure 5. 
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